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(g) Magnetic resonance angiography method and apparatus employing an integration projection. 

(57) Magnetic Resonance Angiography (MRA) 
method and apparatus for determining an 
angiogram, in the form of a two-dimensional 
array of pixel intensities, from time-of-flight or 
phase contrast MRA volumetric image data 
utilizes a weighted sum, pixel by pixel, of inte- 
gration projection (IP) and maximum intensity 
projection (MIP) 2D pixel arrays formed in the 
same viewing direction. The forming of the 
integration projection pixel array is preceded by 
optional background suppression and optional 
application of a depth cueing intensity weight- 
ing function varying along the viewing direc- 
tion. The background suppression procedure 
contains a non-linear unsharp filtering portion 
followed by an intensity thresholding step. 
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BACKGROUND OF THE INVENTION 

1 . Field of the Invention 

The present invention relates to an angiography 
method and apparatus wherein a two-dimensional ar- 
ray of pixels for display as an angiogram is deter- 
mined by a projection or ray casting with respect to a 
volumetric data set obtained from signal samples of 
radiation exiting a region of a body under examina- 
tion, which volumetric data set exhibits a contrast be- 
tween flowing blood and stationary tissue. In its par- 
ticular aspects, the prese/rt invention relates to an- 
giography employing an integration projection (IP). 

2. Background of the Invention 

In the medical imaging field, a three-dimensional 
image data set or volumetric image may be collected 
from X-ray Computed Tomography (CT), Magnetic 
Resonance Imaging (MRI), Positron Emitted Tomog- 
raphy (PET), Single Photon Emission Computed To- 
mography (SPECT), Nuclear Medicine Tomography 
and Ultrasound modalities. In addition, it is a common 
practice to collect a volumetric image as two-dimen- 
sional image data for each of a stack of relatively thin 
slices. However collected, a volumetric image has in- 
tensity values or voxels centered at lattice points in a 
generally regular 3D grid or array. Various of these 
modalities are useful to produce an angiogram from 
a collected volumetric image data set having intensity 
values exhibiting an enhanced contrast between 
flowing blood and stationary tissue, some by employ- 
ing a preparatory step of injecting a contrast agent 
into the bloodstream. 

Magnetic Resonance Angiography (MRA) is car- 
ried out by operating magnetic resonance imaging 
(MRI) apparatus to collect a three-dimensional data 
set of voxels using an RF and magnetic field gradient 
pulse sequence which produces an enhanced con- 
trast between flowing blood and stationary tissue due 
to flow or velocity effects. Commonly used RF and 
magnetic field gradient sequences for angiography 
involve 2D time-of-flight or 2D or 3D phase-contrast 
methods. In the 2D time-of-flight method, a collection 
of spin resonance signals for multiple parallel slices is 
obtained. The flow sensitive contrast is due to sub- 
stantially saturating stationary spins in the slice from 
which spin resonance signals are collected by rela- 
tively rapid repetition of close to 90° flip angle slice se- 
lective RF excitation pulses so that only the unsatu- 
rated spins in blood flowing into the slice have rela- 
tively strong longitudinal magnetization just prior to 
experiencing the excitation pulses. This induces high 
intensity spin resonance signals from the inflowing 
blood, which intensity increases with the amount of in- 
flow velocity component normal to the slice. A three- 
dimensional data set of voxel intensities is computed 


as the amplitudes obtained by two-dimensional Four- 
ier transformation for each slice of samples of the spin 
resonance signals received during a read gradient for 
sequences repeated with different phase encoding 
5 gradient integrals. In the phase contrast method, bi- 
polar gradient pulses are used to produce spin reso- 
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gradient direction. In the 3D method, a three- 
dimensional data set of spin resonance samples is 
10 collected. From these samples a three-dimensional 
data set of voxel intensities is obtained from phase 
values produced by a three-dimensional Fourier 
transformation. However obtained, the initial three- 
dimensional data set can be expanded by interpolat- 
es ing additional intervening layers or slices of voxels to 
. produce a three-dimensional data set to be rendered 
for viewing purposes into a two-dimensional data set 
of pixels. 

Display algorithms for processing three- 

20 dimensional Magnetic Resonance Angiography 
(MRA) data into two-dimensional form are known 
from H. Cline et al, "Volume Rendering and Connec- 
tivity Algorithms for MR Angiography" Magn. Res. 
Med. 18, pp. 384-394 (1991). A typical rendering 

25 method involves the formation of a projection image 
in a projection or viewing direction. Projection meth- 
ods generally involve projecting or casting parallel 
rays through the three-dimensional data set in a view- 
ing direction, the rays being in one-to-one association 

30 with the pixels of the projection image, and from the 
intensities of the voxels intercepted by or interpolated 
or resampled along the respective rays, determining 
the values of the associated pixels. Further, the two^ 
dimensional pixel array determined by the projection 

35 rendering method can be expanded in size prior to 
viewing by interpolating additional intervening rows 
and/or columns of pixels. It is also useful to produce 
projection images from different angles of view cor- 
responding to slowly rotating the viewing direction 

40 and to display these images sequentially In a cine 
loop. Such a cine helps the person reading the angio- 
gram resolve vessel overlap and depth ambiguities. 

The most widely used projection method is max- 
imum intensity projection (MIP). With this method, 

45 which is computationally fast, the maximum intensi- 
ties of voxels along or interpolated along ("interpolat- 
ed" is herein meant to include analogous operations 
such as resampling) the respective rays are taken as 1 
the intensities of the associated pixels. MIP images 

so discard much information that is present in the origi- 
nal volumetric imag e. A single MIP image does not 
provide any depth cues or vessel overlap informatlonT _ 
Consequently, it is difficult to resolve depth ambigu- t^u^c 
ities and vessel overlap even from a cine. Further, 

55 there is a loss of signal at vessel edges which causes 
apparent vessel narrowing and loss of vessel thick- 
ness information. Another deficiency of this method is 
that tortuous and looping vessels are difficult to dis- 
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tinguish from vascular aneurysms. 

One alternative projection method is the so-cal- 
led integration projection (IP) method in which the 
sums of the intensities of voxels along (or interpolated 
along) the respective rays are taken as the values of 5 
the associated pixels. Fast implementations can be 
achieved by using the relationship that a 2D integra- 
tion projection image of a 3D data set in a projection 
direction is proportional to the inverse 2D Fourier 
transform of the 2D data set intercepted by (or inter- 10 
polated along) a plane through the origin of the Four- 
ier transform of the 3D data set whose normal is in the 
projection direction. The IP method preserves vessel 
thickness information because pixel intensity increas- 
es with increased vessel thickness in the projection is 
direction. Unfortunately, there are also large contri- 
butions from integrated background voxels which ob- 
scure important vessel detail in the angiogram image. 
A known variation of the IP method is to sum only 
those voxels which have intensities above a given in- 20 
tensity threshold set by the operator. If the threshold 
is set too low, the contribution of background voxels 
still obscure vessel detail and result in a generally noi- 
sy image. Conversely, if the threshold is set too high, 
there is a loss of small vessel and edge information. 25 

P.J. Keller et al., "SIMP: An Integrative Combina- 
tion with MIP'\ Proceedings of the Annual Meeting of 
the society of Magnetic resonance in Medicine, 
(1991) p.201 discloses a combination of IP and MIP 
methods. A scaled integration projection is used for 30 
each ray intercepting at least two voxels that have in- 
tensities above a user-defined threshold. For rays not 
intercepting at least two such voxels, a maximum in- 
tensity projection is used. Too low a threshold results 
in an angiogram primarily having the characteristics 35 
of an IP image, while too high a threshold results in 
an angiogram primarily having the characteristics of 
a MIP image. 

F.R. Korosec et al., "A Data Adaptive Reprojec- 
tion Technique for MR Angiography", Proceedings of 40 
the Annual Meeting of the Society of Magnetic Reso- 
nance in Medicine", (1991), p.820, discloses another 
method that combines integration and maximum in- 
tensity projections. In this method, the original volu- 
metric image data is compared with a first intensity 45 
threshold to def ine a vessel mask which is then blur- 
red and converted to binary form. The blurred binary 
mask is applied to the original volumetric image data 
and thereafter a second, presumably higher, intensity 
threshold is applied to form processed volumetric im- so 
age data. The final image is a weighted average of an 
integration projection image of the processed volu- 
metric image date, after roll-off by a non-linear inten- 
sity transformation, and a maximum intensity projec- 
tion of the original volumetric image data. This meth- . 55 
od appears critically dependant on the values set for 
the thresholds and further, due to the absence of 
depth cues, ambiguities remain as to the relative pos- 


itions of overlapping vessels, which must be resolved 
from sequential projections at incremental angles. 
Further, setting the thresholds sufficiently high to 
suppress background may eliminate vessel detail in 
the integration projection contribution to the final im- 
age and such detail may be too faint, without integra- 
tion along the direction of view, to be apparent in the 
maximum intensity projection contribution to the final 
image. 

J. Listerud et al., "TRAP: Traced Ray by Array 
Processor, Proceedings of the Annual Meeting of the 
Society of Magnetic Resonance in Medicine, (1991) 
p.757, describes pseudo-surface rendering utilizing 
dynamic range compression and artificial lighting fol- 
lowed by maximum intensity projection. Artificial 
lighting of a range image makes structures closer to 
the viewer brighter than and appearing superimposed 
in front of those which are more distant. However, the 
maximum intensity projection of the range image de- 
stroys much of the depth information produced by the 
artificial lighting as well as suffering from the already 
described deficiencies of MIP. 

SUMMARY OF THE INVENTION 

It is an object of the present invention to provide 
an angiography method and apparatus for producing 
an angiogram containing details important for inter- 
pretation including vessel overlap and depth cue in- 
formation, which are not obscured by background 
noise. 

It is also an object of the present invention to pro- 
vide an angiography method and apparatus which 
employs depth cueing followed by an integration pro- 
jection. 

It is yet another object of the present invention to 
provide an angiography method and apparatus which 
employs background suppression followed by an in- 
tegration projection, which background suppression 
preserves the edges of vessel details. 

It is still a further object of the present invention 
to provide an angiography method and apparatus 
combining integration and maximum intensity projec- 
tions which includes depth cueing and/or background 
suppression prior to the integration projection. 

Briefly, the aforementioned and other objects of 
the present invention are satisfied by trie provision of 
an angiography method and apparatus for determin- 
ing a two-dimensional array of pixel intensities, in 
which radiation, for example MR spin resonance RF 
signals, is induced to exit from volume elements of a 
region of a body under examination, the radiation hav- 
ing a characteristic, for example - ampl itude or phase, 
influenced by any blood within said volume elements, 
said flowing blood being in cardiovascular structure 
including blood vessels, the volume elements having 
center to center spacings in respective three linear in- 
dependent, for example, mutually orthogonal direo 
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tions. The radiation exiting from the volume elements 
of the region is received and sampled to form a col- 
lection of signal samples. This collection of signal 
samples is converted into an initial three-dimensional 
array of voxel intensities, each dimension of which 
corresponds to a different one of the three mutually 
orthogonal directions. Afurtherthree-dimensional ar- 
ray, herein referred to as a pre-depth-cued array, is 
derived therefrom by performing desired preliminary 
operations, if any, including but not limited to expan- 
sion of the size of the array through interpolation of 
additional intervening slices or layers and/or perform- 
ing a background suppression operation thereon. The 
voxel intensities of this pre-depth-cued three- 
dimensional array are modified in accordance with a 
depth cueing function applied in a given direction to 
form a three-dimensional depth-cued array of voxel 
intensities. Then, the intensities of voxels intercepted 
by or interpolated along respective parallel rays pro- 
jected through the depth-cued array of computed vox- 
el intensities in the given direction, to form respective 
intensity sums for the respective rays, each of the 
rays being associated with a different one of the pix- 
els of the two-dimensional array to be determined. 
Lastly, the values of the respective pixels of the two- 
dimensional array are determined as a function of the 
respective sums for the rays associated with the re- 
spective pixels. 

In accordance with a further aspect of the present 
invention, in order to produce an angiogram in which 
integrated background voxels are inhibited from ob- 
scuring vessel detail while also preserving the edges 
of the vessels, a background suppression procedure 
is applied to a pre- background-suppressed three- 
dimensional array derived from the initial three- 
dimensional array by desired preliminary operations, 
if any, including but not limited to the aforementioned 
array size expansion by interpolation in order to pro- 
duce background suppressed three-dimensional ar- 
ray. The background suppression procedure compris- 
es a non-linear unsharp filtering portion followed by 
an intensity thresholding step. The non-linear un- 
sharp filtering portion comprises comparing the inten- 
sities of the voxels of the pre-backg round-sup- 
pressed array with a first intensity value, namely, a 
maximum background intensity value, to segment the 
voxels into flow voxels and background voxels, form- 
ing a clipped array of voxels by replacing the flow vox- 
els of the pre-background suppressed array with a 
second intensity value less than the first intensity val- 
ue, namely, the mean of the intensities of the back- 
ground voxels, applying a low-pass spatial filter to the 
clipped array to form a smoothed clipped array, and 
subtracting the smoothed clipped array from the pre- 
backg round-suppressed array, voxel by voxel. The in- 
tensity thresholding step forms a pre-depth-cued or 
background suppressed array. This step utilizes, as 
a threshold, a third intensity value computed as the 


first intensity value minus the second intensity value. 

A further aspect of the present invention involves 
the determination of the respective maximum inten- 
sities of voxels intercepted by or interpolated along 

5 the respective parallel rays associated with the re- 
spective pixels of the two-dimensional array and di- 
rected in the samfi direction as the depth cueing func- 
tion, as projected through a pre-MIP three- 
dimensional array of voxel intensities. The three- 

10 dimensional pre-MIP array is derived from the initial 
three-dimensional array by preliminary operations, if 
any, including but not limited to the aforementioned 
expansion of array size by interpolation. Further, the 
aforementioned determination of the values of the re- 

15 spective pixels of the two-dimensional array is by 
forming a weighted combination of the respective 
maximum intensities and the respective sums for the 
rays associated with the respective pixels. 

The integration projection after background sup- 

20 pression and depth cueing presents vessel details in- 
cluding vessel overlap in a manner which is readily in- 
ter pretable. The combination therewith of the maxi- 
mum intensity projection is useful to introduce details 
that might be overly suppressed in forming the inte- 

25 gration projection. 

Other objects, features and advantages of the 
present invention will become apparent upon review 
of the following detailed description, when taken in 
conjunction with the appended drawing, wherein: 

30 

BRIEF DESCRIPTION OF THE DRAWING 

Figure 1 is a simplified schematic of magnetic 
resonance imaging apparatus, including an im- 
35 age processor, usable for implementing the pres- 

ent invention; 

Figure 2 is an illustration of the production of an 
angiogram by projection; 

Figure 3 is a flow chart of the operations, includ- 
40 ing background suppression and depth cueing 

procedures, performed by the image processor 
of Figure 1; 

Figure 4 is a flow chart of the background sup- 
pression procedure of Figure 3; 
45 Figure 5 is a typical histogram of the input array 

to the background suppression procedure of Fig- 
ure 3; and 

Figure 6 is a graph of a depth cueing function 
used in the depth cueing procedure of Figure 3. 

50 

DETAILED DESCRIPTION OF THE PREFERRED 
EMBODIMENT 

Figure 1 of the drawing depicts a simplified sche- 
55 matic of a magnetic resonance imaging (MRI) appa- 
ratus 1 0, which is given for the purpose of orientation. 
Therein, a main field coil 12 establishes a steady unh 
form longitudinally directed magnetic field in an elon- 
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gated measurement space 14. The approximate 
range of field strengths in present MRI equipment is 
in the range of 0.1 to 4.0 Tesla. Main field coil 12 may 
be resistive or superconducting, in the latter case 
power supply 16 is used solely for initiating the nec- 5 
essary coil current Also acting upon measurement 
space 14 are a system of three independently control- 
lable gradient field coils 18 producing longitudinally 
directed fields with strengths that vary linearly in re- 
spective three orthogonal directions, in proportion to 10 
the respective coil currents driven by power supply 
20. An RF coil 22 receives RF pulses formed in RF 
transmitter section 24 via Transmit/Receive switch 26 
and transmits these pulses into measurement space 
14. As is well known, these RF pulses are at the Lar- 15 
mor frequency determined by the strength of the main 
field and the nuclei of interest (usually water protons). 
RF coil 22 receives the resultant spin resonance sig- 
nals from the measurement space and via switch 26 
feeds these signals to receiver section 28 where they 20 
are detected and then sampled and digitized by an 
A/D converter 30. A sequence controller 31 controls 
gradient power supply 20, Transmit/Receive switch 
26, RF transmitter section 24 and RF receiver section 
28 so as to subject measurement space 14 to a se- 25 
quence of sub-experiments, each including at least a 
read or excitation RF pulse for converting some or all 
longitudinal magnetization to transverse magnetiza- 
tion, and gradient pulses functioning as one or more 
phase encoding gradients (the time integral of which 30 
is stepped from each sub-experiment to the next) fol- 
lowed by a read gradient during which spin resonance 
signals are sampled by the A/D converter 30. 

Typically, for angiography using a time-of-f light 
procedure, the sequence is configured to induce spin 35 
resonance signals of relatively high amplitude from 
voxels through which blood is flowing than from vox- 
els containing only stationary tissue. A so-called 2D 
"multislice" collection is commonly utilized in which a 
slice-selective sequence is applied sequentially to a 40 
stack of a plurality of (for example, five) equally 
spaced apart parallel slices, of known thickness. In a 
slice-selective sequence; at least the RF read pulse 
is applied simultaneously with a slice-selection gra- 
dient. The result of the multislice collection is a two- 45 
dimensional array (containing, for example, 256 X 
256 points in a spatial frequency space known as K- 
space) of signal samples for each slice, so that the 
multi-slice collection is effectively a three- 
dimensional array of signal samples which are sup- so 
plied by A/D converter 30 to digital memory 32. The 
signal sample data for each slice is separately sub- 
jected to a two-dimensional Fourier Transform in 
block 34, which may constitute an array processor, to 
produce a two-dimensional array of pixel intensities 55, 
for each slice. As a result, the collection of pixels for 
multiple parallel slices supplied by block 34 to digital 
memory 36 can be considered a three-dimensional 


array of voxel intensities (of, for example, 5 X 256 X 
256 points). - 

In the case of the use of either a 3D time-of-f light 
method or 3D phase-sensitive method, a three- 
dimensional array of sample points in k-space is read 
into memory 32 with usually one of the dimensions, 
which is referred to herein as vertical, being quan- 
tized by less points than the other dimensions. This 
three-dimensional k-space array is then subjected in 
block 34 to a three-dimensional Fourier Transform to 
produce the three-dimensional array of voxel intensi- 
ties in memory 36. 

This three-dimensional array of voxel intensities 
is operated upon by an image processor 38 to supply 
to digital memory 40 a two-dimensional array of pixel 
intensities which are displayed on a display device 42, 
such as a CRT. Image processor, in accordance with 
the present invention, includes program modules for 
various functions including background suppression 
38a, interpolation 38b, projection 38c and depth cue- 
ing 38d, which are pertinent to angiography. 

Figure 2 illustrates application of projection mod- 
ule 38c to produce a two-dimensional array of pixels 
for display as an angiogram 44 which is obtained by 
projecting parallel rays in a viewing direction 46 
through a three-dimensional array 48 of voxels which 
are pixels of five stacked slices 50a through 50e, 
each ray producing an associated pixel of the angio- 
gram. For purposes of illustration, a three- 
dimensional intensity pattern corresponding to flow- 
ing blood 52 within a vessel is positioned intersecting 
three-dimensional array 48. As a result of the projec- 
tion procedure, the voxels making up the disc shaped 
volumes 54a, 54b of vessel branches 52a, 52b within 
slice 50a, project respectively to pixels making up the 
areas 55a, 55b of angiogram 44. The most common 
projection technique for MR angiography is maximum 
intensity projection (MIP). In MIP, the intensity of each 
pixel in the angiogram is taken as the maximum inten- 
sity of the voxels through which the associated ray 
passes. A less common projection technique is inte- 
gration projection (IP), wherein, the intensity of each 
pixel of the angiogram is taken as the sum of the in- 
tensities of the voxels through which the associated 
ray passes. As, is apparent in the example given, in 
the absence of interpolation either before or after pro- 
jection, the angiogram produced in a horizontal view- 
ing direction of the illustrated stack of slices 50a 
through 50e would be only five pixels high. 

The operation of image processor 38 is indicated 
by the flow chart of Figure 3. The volumetric data set 
56 stored in memory 36 is the input thereto and the 
angiogram 57 stored in memory 40 is the output 
therefrom. First it is determined at step 58 whether 
the number of voxels in the volumetric data set is to 
be expanded. This may be determined by operator in- 
put, but preferably a default condition is provided in 
which whether or not to interpolate, and the number 
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and orientation of intervening slices or layers to be in- 
terpolated, if any, is given by default based upon the 
spatial configuration of the volumetric data set 56 pre- 
sented. Interpolation procedure 60 is performed by 
module 38b only if the determination at step 58 is 5 
YES. As noted, if only a small number of slices are ac- 
quired, it would be desirable to improve apparent re- 
solution by interpolating additional voxels at least in 
the slice direction. The preferred technique for inter- 
polation of each slice intervening originally adjoining 10 
slices (or layers intervening originally adjoining lay- 
ers) is setforth in EP-A-0 571 026 entitled "Direction- 
al Interpolation For Magnetic Resonance Angiogra- 
phy" having a priority date of May 21, 1992, which is 
herein incorporated by reference. In the directional in- is 
terpolation procedure disclosed and claimed therein, 
there is a determination for each voxel in one of the 
originally adjoining slices whether a voxel exists in the 
other of the originally adjoining slices which corre- 
sponds thereto with respect to position relative to lo- 20 
cal intensity patterns. A local direction of interpolation 
directed between these voxels is used if correspon- 
dence is found. Otherwise, the conventional local di- 
rection of interpolation normal to the slice is used. 
Such directional interpolation eliminates the staircase 25 
artifact which may otherwise appear along edges of 
obliquely directed vessels. 

Next, it is determined at step 62 whether back- 
ground is too high, necessitating background sup- 
pression procedure 64 to be performed by module 30 
38a. This procedure is detailed in Figure 4 and in Ta- 
ble I and will be seen to comprise a nonlinear unsharp 
filtering procedure portion comprising steps 64a-64f 
followed by an intensity thresholding step 64g. The 
determination of whether to perform background sup- 35 
pression also may be by operator input. In step 64b, 
the voxels of the input array are roughly segmented 
into flow voxels and background voxels and a mean 
intensity of the background voxels is determined. This 
segmentation is accomplished by comparing voxel in- 40 
tensities with a maximum background intensity, f jow 
voxels having intensities greater than the maximum 
background intensity and background voxels having 
intensities less than or equal to the maximum back- 
ground voxel intensity. 45 

While the maximum background intensity could 
be a value which is input by an operator, preferably 
the maximum background intensity is automatically 
determined in step 64a for an operator designated Re- 
gion of Interest (ROI) containing blood vessel struc- so 
ture. Atypical histogram of such a region is shown in 
Figure 5 as having a high frequency of occurrence or 
number of voxels peak at a low intensity approximate- 
ly at a mean background intensity l^mean, and smaller 
frequency of occurrence peaks at intensities greater 55 
than a maximum background intensity l^max. corre- 
sponding to flow. If l moan and l^ are the mean and 
standard deviation, respectively, of the intensities of 
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all the voxels in the ROI, then the maximum back- 
ground intensity is determined in step 64a, as follows: 

lbg_max = 'mean + C • l^ 

where C is a constant. For time-of flight images, it 
was found empirically that a value of C = 2.0 is a good 
compromise between suppressing the background 
without losing the small vessels. .After the afcrerners 
tioned segmentation, the mean intensity of the back- 
ground, Ibgumeam is determined in step 64c as the aver- 
age of all voxel intensities within the ROI having an 
intensity less than or equal to l^max- In step 64d, a 
clipped array is formed by setting the intensities of all 
voxels corresponding to flow in the array input to pro- 
cedure 64 (those having intensities greater than the 
maximum background intensity lbQ_ max ) to the mean 
background intensity l^m^n- Then, in step 64e a low 
pass filtering kernel is applied to the clipped array to 
form a smoothed clipped array. A boxcar averaging 
filter which is 3 voxels wide in each dimension filtered 
is a suitable kernel. Thereafter, in step 64f, the 
smoothed clipped array, is subtracted from the input 
array to procedure 64. The aforementioned non-linear 
unsharp filtering portion effects a high pass filtering 
operation that does not blur vessel edges. It generally 
reduces the range of background intensities to zero 
through Ibgumax-lbajnean- Then, in intensity thresholding 
step 64g, voxels having an intensity less than or equal 

to 'thresh = ^bg_ max-lba_mean are Set to Zero. 

In step 66 it is determined by operator input 
whether depth cues are to be provided. If the depth 
cueing option is selected, then in step 68 a scale func- 
tion having an intensity multiplying scale factor that 
monotonically increases over distance in the viewing 
direction from back to front depth reference locations 
fixed with respect to the volumetric image is used to 
introduce a contrast whereby closer objects appear 
brighter than more distant ones. Although various 
monotonic functions are useful, that currently consid- 
ered preferable is graphed in Figure 6 and will be 
seen to be a linearly increasing function. 

Referring to Figure 6, assume for the purpose of 
illustration that the viewing direction is along an axis 
Y. Front and back depth cue reference planes are in- 
dicated at 

Y = Y f and Y = Y b , respectively. Then the scale factor, 
S, for multiplying the intensity of a voxel that is be- 
tween the two planes at location Y is given by: 

S = S b + (Y - Y b )#(Sf - S b )/(Y f - Y b ) 
where the scale factors for voxels at the front and 
back depth cue reference planes are S = Sf and S = 
S b , respectively. For locations Y behind Y b , S = S b , 
and for locations Y in front of Yf, S = Sf. 

Following the optional depth cueing procedure, 
an integration projection is performed in step 70 by 
assigning to the pixels associated with parallel rays 
projected through the volumetric image in a viewing 
direction (generally illustrated in Figure 2), the sum of 
the intensities of the voxels intercepted bytherespec- 
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tive rays. It should be understood that when, in a pro- 
jection procedure, the parallel rays do not pass pre- 
cisely through the centers of the voxels, it is conven- 
tional to interpolate voxel intensities at points along 
the rays, or equivalently, to resample (by an interpo- 
lation process) the volumetric image to determine 
new voxel centers along the rays. A tri-linear interpo- 
lation of voxel intensities along projection rays is well 
known and preferred. 

The resultant two-dimensional image is then 
scaled in intensity to lie in the range zero to l newma x. 
by using the following linear formula: 

'new = 'newmax # Ootd " 'thresh)^ old max .." 'thresh) 

where: 

loid = Intensity of a pixel having a non-zero Intensity 
prior to scaling. 

I new - Intensity of the pixel after scaling, 
lnewmax = A maximum allowed intensity value after 
scaling (usually 255 for byte scaled images), 
^otdmax = Maximum intensity of the pixels prior to scal- 
ing. 

'thresh " The intensity threshold used in step 64g. 
Those pixels having a zero intensity prior to scaling 
are not treated by the scaling operation and therefore 
retain a zero intensity after scaling. 

The resultant scaled integration projection (IP) 
image shows vessel overlap and thickness. When the 
background suppression procedure 64 has been em- 
ployed, the IP image has lower background levels 
than a conventional maximum intensity projection 
(MIP) image. Further, when the depth cueing step 68 
has been used, depths cues enable visualization of 
the relative positions of overlapping vessels. The IP 
image compares unfavourably to an MIP image in re- 
spect to small vessels. These appear narrow and 
dark. It was found that the deficiencies of the resultant 
IP image can be compensated by the introduction of 
a further pixel intensity component due to maximum 
intensity projection. 

Accordingly, in step 74, a maximum intensity pro- 
jection is performed with respect to the pre-back- 
ground -suppression volumetric image (the volumet- 
ric image input to background suppression selection 
step 62) in the same direction and with the same con- 
figuration of parallel projection rays as used in the in- 
tegration projection (IP) step 70. As is conventional, 
the intensities of the pixels associated with the re- 
spective projection rays are determined as the maxi- 
ma of the voxel intensities intercepted by or interpo- 
lated along the respective rays. Then in step 76, a 
weighted sum, pixel by pixel, of the resultant MIP im- 
age and the scaled IP image is formed in accordance 
with the following formula: 

Uum = Omlp + D#l new )/(1 + D) 

where: 

l m ip = Intensity of a pixel in the MIP image, 
'new = Intensity of the same pixel in the scaled IP im- 
age. 


D = A constant preferably in the range of 1 to 4. 

Prior to forming the weighted sum in step 76, in 
step 72 the intensities, l m{p , of the MIP image are also 
scaled to be in the range zero to lnewmax- 

5 After the formation of the weighted sum two-di- 

mensional image, it is determined in step 78, for ex- 
ample, from a number of lines or pixels to be spanned 
fc by the ROI in display 42, whether to expand the num- 
ber of pixels in the weighted sum image. If an expand- 
ed ed number of pixels is indicated, then this is done by 
the 2D interpolation step 80, to determine further 
lines of pixels between existing lines. Preferably, this 
step also uses the directional interpolation techni- 
ques in my prior patent application which was herein 

15 previously referred to and incorporated by reference, 
in connection with the 3D interpolation step 60. 

The final angiogram image 57 has the small ves- 
sel detail of the MIP image, and the vessel crossing, 
thickness, and depth information of the depth-cued 

20 integration projection image. The background sup- 
pression procedure 64 has been found useful for 
time-of-flight volumetric image data, where back- 
ground values are high. Phase contrast volumetric 
image data, however, has high vessel -to- background 

25 contrast and may not require, background suppres- 
sion. The depth cueing step has been found very use- 
ful for visualizing tortuous cerebral arteries. 

The present invention has been found to provide 
more information in each projection and to therefore 

30 require fewer projections from different directions of 
view to enable interpretation of the angiogram. Stud- 
ies for which the method appears promising is the cir- 
cle of wiilis, sinus, carotid bifurcation, carotid siphon, 
pulmonary and peripheral vessels. 

35 It should now be apparent that the objects of the 

present invention have been satisfied, and further 
that while the present invention has been described 
in detail with respect to MR angiography, many prin- 
ciples thereof are applicable to other medical imaging 

40 modalities including but not limited to CT and nuclear 
medicine tomography. Further, numerous modifica- 
tions in the details given are possible within the in- 
tended spirit and scope of the invention. 

45 

Claims 

1. An angiography method for determining a two-di- 
mensional array of pixel intensities, said method 

so comprising: 

inducing radiation to exit from volume ele- 
ments of a region of a body under examination, 
said radiation having a characteristic influenced 
by any flowing blood within said volume ele- 

55 ments, said blood being in cardiovascular struc- 

ture including blood vessels, said volume ele- 
ments having center to center spacings in re- 
spective three linear independent directions; 
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receiving and sampling the radiation exit- 
ing from the volume elements of said region to 
form a collection of signal samples; 

converting said collection of signal sam- 
ples into an initial three-dimensional array of vox- 5 
el intensities, each dimension of said array corre- 
sponding to a different one of said three linear in- 
dependent directions; 

modifying the voxel intensities of a pre- - 
depth-cued three-dimensional array of voxel in- 10 
tensities derived from said initial array in accor- 
dance with a depth cueing function applied in a 
given direction to form a three-dimensional 
depth-cued array of voxel intensities; 

summing the intensities of voxels inter- 15 
cepted by or interpolated along respective paral- 
lel rays projected in said given direction through 
said depth-cued array of computed voxel intensi- 
ties to form respective intensity sums for the re- 
spective rays, each of said rays being associated 20 
with a different one of the pixels of the two-di- 
mensional array to be determined; and 

determining the values of the respective 
pixels of said two-dimensional array as a function 
of the respective sums for the rays associated 25 
with the respective pixels. 

A method as claimed in Claim 1 , further compris- 
ing determining respective maximum intensities 
of voxels intercepted by or interpolated along said 30 
respective parallel rays associated with the re- 
spective pixels and directed in said given direc- 
tion, as projected through a pre-MIP three- 
dimensional array of voxel intensities derived 
from said initial three-dimensional array, and 35 
wherein said determining of the values of the re- 
spective pixels of said two-dimensional array is 
by forming a weighted combination of said re- 
spective maximum intensities and said respec- 
tive sums for the rays associated with said re- 40 
spective pixels. 

A method as claimed in Claim 1 or 2, wherein a 
pre-background-suppressed three-dimensional 
array is derived from said initial three- 45 
dimensional array and said pre-depth-cued 
three-dimensional array is derived from said pre- 
background-suppressed array by the following 
operations: 

comparing the intensities of the voxels of 50 
said pre-background-suppressed array with a 
first intensity value to segment said voxels into 
flow voxels and background voxels; 

forming a clipped array of voxels by replac- 
ing the flow voxels of said pre- background sup- 55 
pressed array with a second intensity value less 
than the first intensity value; 

applying a low-pass spatial filter to said 


clipped array to form a smoothed clipped array; 

subtracting said smoothed clipped array 
from said pre-background suppressed array, 
voxel by voxel to form a pre-thresholded array; 
and 

intensity thresholding said pre-threshold- 
ed array with respect to a third intensity value to 
form said pre-depth-cued array. 

4. An angiography method for determining a two-di- 
mensional array of pixel intensities, said method 
comprising: 

inducing radiation to exit from volume ele- 
ments of a region of a body under examination, 
said radiation having a characteristic influenced 
by any blood within said volume elements, said 
blood being in cardiovascular structure including 
blood vessels, said volume elements having cen- 
ter to center spacings in respective three linear in- 
dependent directions; 

receiving and sampling the radiation exit- 
ing from the volume elements of said region to 
form a collection of signal samples; 

converting said collection of signal sam- 
ples into an initial three-dimensional array of vox- 
el intensities, each dimension of said array corre- 
sponding to a different one of said three linear in- 
dependent directions; 

deriving a background suppressed three 
dimensional array by the following procedure 
from a pre-background-suppressed three- 
dimensional array derived from said initial three- 
dimensional array: 

a) comparing the intensities of the voxels of 
said pre-background-suppressed array with a 
first intensity value to segment said voxels 
into flow voxels and background voxels; 

b) forming a clipped array of voxels by replac- 
ing the flow voxels of said pre-background 
suppressed array with a second intensity var- 
ue less than the first intensity value; 

c) applying a low-pass spatial filter to said 
clipped array to form a smoothed clipped ar- 
ray; 

d) subtracting said smoothed clipped array 
from said pre-background suppressed array, 
voxel by voxel to form a pre-thresholded ar- 
ray; and 

e) intensity thresholding said pre-thresholded 
array with respect to a third intensity value to 
form said background suppressed array; 

summing the intensities of voxels inter- 
cepted by or interpolated along respective paral- 
lel rays projected in a given direction through a 
pre-projection three dimensional array derived 
from said background suppressed array, to form 
respective intensity sums for the respective rays, 
each of said rays being associated with a differ- 
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entone of the pixels of the two-dimensional array 
to be determined; and 

< determining the values of the respective 
pixels of said two-dimensional array as a function 
of the respective sums for the rays associated 
with the respective pixels. 

5. A method as claimed in Claim 4, further compris- 
ing determining respective maximum intensities 
of voxels intercepted by or interpolated along said 
respective parallel rays associated with the re- 
spective pixels and directed in said given direc- 
tion, as projected through a pre-MIP three- 
dimensional array of voxel intensities derived 
from said initial three-dimensional array, and 
wherein said determining of the values of the re- 
spective pixels of said two-dimensional array is 
by forming a weighted combination of said re- 
spective maximum intensities and said respec- 
tive sums for the rays associated with said re- 
spective pixels. 

6. A method as claimed in Claim 1, 2, 3, 4 or 5 
wherein said induced and received radiation are 
MR spin resonance RF signals. 

7. An angiography apparatus for determining a two- 
dimensional array of pixel intensities, said appa- 
ratus comprising: 

means for inducing radiation to exit from 
volume elements of a region of a body under ex- 
amination, said radiation having a characteristic 
influenced by any blood within said volume ele- 
ments, said blood being in cardiovascular struc- 
ture including blood vessels, said volume ele- 
ments having center to center spacings in re- 
spective three linear independent directions; 

means for receiving and sampling the ra- 
diation exiting from the volume elements of said 
region to form a collection of signal samples; 

means for converting said collection of sig- 
nal samples into an initial three-dimensional ar- 
ray of voxel intensities, each dimension of said 
array corresponding to a different one of said 
three linear independent directions; and 

image processing means comprising: 

a) means for modifying the voxel intensities of 
a pre-depth-cued three-dimensional array of 
voxel intensities derived from said initial array 
in accordance with a depth cueing function 
applied in a given direction to form a three- 
dimensional depth cued array of voxel inten- 
sities; 

b) means for summing the intensities of voxels 
intercepted by or interpolated onto respective 
parallel rays projected through said depth 
cued array of computed voxel intensities in 
said given direction, to form respective inten- 


sity sums for the respective rays, each of said 
rays being associated with a different one of 
the pixels of the two-dimensional array to be 
determined; and 
5 c) means for determining the values of the re- 

spective pixels of said two-dimensional array 
as a function of the respective sums for the 
rays associated with the respective pixels. 

10 8. An apparatus as claimed in Claim 7, wherein said 
image rendering means further comprises: 

d) means for determining respective max- 
imum intensities of voxels intercepted by or inter- 
polated along said respective parallel rays asse- 
ts ciated with the respective pixels and directed in 
said given direction, as projected through a pre- 
MIP three-dimensional array of voxel intensities 
derived from said initialthree-dimensional array; 
and wherein said means for determining 
20 the values of the respective pixels of said two-di- 
mensional array comprises means for forming a 
weighted combination of said respective maxi- 
mum intensities and said respective sums for the 
rays associated with said respective pixels. 

25 

9. An apparatus as claimed in Claim 7 or 8, wherein 
a pre-background-suppressed three-dimensional 
array is derived from said initial three-dimensional 
array and wherein said image processing means 

30 further comprises background suppression 

means for converting said pre-background sup- 
pressed array to said pre-depth-cued array, said 
background suppression means comprising 
means for 

35 comparing the intensities of the voxels of 

said pre-background-suppressed array with a 
first intensity value to segment said voxels into 
flow voxels and background voxels; 

forming a clipped array of voxels by replao 

40 ing the flow voxels of said pre-background sup- 

pressed array with a second intensity value less 
than the first intensity value; 

applying a low-pass spatial filter to said 
clipped array to form a smoothed clipped array; 

45 subtracting said smoothed clipped array 

from said pre-background suppressed array, 
voxel by voxel to form a pre-thresholded array; 
and 

intensity thresholding said pre-threshold- 
50 ed array with respect to a third intensity value to 

form said pre-depth-cued array. 

10. An apparatus as claimed in Claim 7, 8 or 9, 
wherein said induced and received radiation are 

55 MR spin resonance RF signals. 
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